perm filename BBESCO.F4[1,MUS] blob sn#078092 filedate 1973-12-18 generic text, type T, neo UTF8
00100	C This is taken from the "Handbook of Mathematical Functions"
00200	C by Abramowitz and Stegun, Dover press S1272, 1965.
00300	C Chapter 9 is on Bessel Functions. The polynomial
00400	C approximations are found in section 9.4, pages 369
00500	C and 370. The recursion formula for generating higher
00600	C orders from J0 and J1 is found at the bottom of
00700	C table 9.1, page 390
00800	
00900		SUBROUTINE BESCOEF(MI,X,IX,K)
01000		REAL MI,J0,J1
01100		COMMON FREQ(3,0/50,100),FUNC(100),AMP(100),II(1),IJJ(4000)
01200		IF(IX.GT.1)GO TO 30
01300		IF(MI.GE.3.0)GO TO 10
01400		xs=(MI/3)**2
01500		j0=1.0+xs*(-2.2499997+xs*(1.2656208+xs*(-0.3163866+
01600		1xs*(0.0444479+xs*(-0.0039444+xs*0.00021)))))
01700		j1=MI*(.5+xs*(-0.56249985+xs*(0.21093573+
01800		1xs*(-0.03954289+xs*(0.00443319+xs*(-0.00031761+
01900		2xs*0.00001109))))))
02000		GO TO 20
02100	10	xs=3.0/MI
02200		j0=(1.0/sqrt(MI))*(0.79788456+xs*(-0.00000077+
02300		1xs*(-0.0055274+xs*(-0.00009512+xs*(0.00137237+
02400		2xs*(-0.00072805+xs*0.00014476))))))*
02500		3cos(MI-0.78539816+xs*(-0.04166397+
02600		4xs*(-0.00003954+xs*(0.00262573+
02700		5xs*(-0.00054125+xs*(-0.00029333+
02800		6xs*0.00013558))))))
02900		j1=(1.0/sqrt(MI))*(0.79788456+xs*(0.00000156+
03000		1xs*(0.01659667+
03100		1xs*(0.00017105+xs*(-0.00249511+xs*(0.00113653-
03200		2xs*0.00020033))))))*
03300	 	3cos(MI-2.35619449+xs*(0.12499612+xs*(0.0000565+
03400		4xs*(-0.00637879+xs*(0.00074348+xs*(0.00079824-
03500		5xs*0.00029166))))))
03600	20	FREQ(2,0,K)=J0
03700		FREQ(2,1,K)=J1
03800		FREQ(2,2,K)=J1
03900		RETURN
04000	30	Y=IX
04100		IF(MI.LE.0.0000001)GO TO 50
04200		IJ=IX*2		
04300		X=((2.0*(Y-1.0))/MI)*FREQ(2,IJ-2,K)-FREQ(2,IJ-4,K)
04400	40	RETURN
04500	50	X=0
04600		RETURN
04700		END